Gravitational radiation reaction and inspiral waveforms in the adiabatic limit 
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We describe progress evolving an important limit of binary orbits in general relativity, that of a stellar mass 
compact object gradually spiraling into a much larger, massive black hole. These systems are of great interest 
for gravitational wave observations. We have developed tools to compute for the first time the radiated fluxes of 
energy and angular momentum, as well as instantaneous snapshot waveforms, for generic geodesic orbits. For 
special classes of orbits, we compute the orbital evolution and waveforms for the complete inspiral by imposing 
global conservation of energy and angular momentum. For fully generic orbits, inspirals and waveforms can 
be obtained by augmenting our approach with a prescription for the self force in the adiabatic limit derived by 
Mino. The resulting waveforms should be sufficiently accurate to be used in future gravitational-wave searches. 
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The late dynamics of a merging compact binary remains 
one of the greatest challenges of general relativity (GR). GR 
doesn't have a "two body" problem so much as it has a "one 
spacetime" problem: one must use the Einstein field equations 
to find the dynamical spacetime describing a multibody sys- 
tem. Although numerical relativity has made great progress in 
recent years (e.g., 1 11] and references therein), most astrophys- 
ically relevant progress has come from identifying a small pa- 
rameter which defines a perturbative expansion. One such ap- 
proach is post-Newtonian (PN) theory (e.g., 01 and references 
therein) — essentially, an expansion in interaction potential 
GM/rc 2 . PN theory works very well when the bodies' sep- 
aration r is large. As the bodies come close, the expansion 
must be iterated to high order (though it may be possible to 
modify the expansion to improve its convergence 1 3]). 

Strong field binaries can be modeled very accurately if 
one member is far more massive than the other. The space- 
time is then that of the larger body plus a perturbation due 
to the smaller body, with the mass ratio acting as an expan- 
sion parameter. This limit is not just of formal interest: bi- 
naries consisting of "small" compact bodies (white dwarfs, 
neutron stars, or black holes with mass fi ~ 1 — 100 Mq) 
captured onto highly eccentric orbits of massive black holes 
(M ~ 10 5 — 1O 7 M ) are important targets for space-based 
gravitational-wave (GW) antennae, especially the LISA |@] 
mission. Such captures are estimated to occur with a rate den- 
sity ~ 10 ±1 Gpc~ 3 yr _1 . Convolving with LISA'S planned 
sensitivity, one expects to measure dozens to thousands of 
events per year |5]. Precisely measuring GW phase as the 
smaller body spirals into the black hole (driven by GW back- 
reaction) and fitting to detailed models will determine sys- 
tem parameters with extraordinary accuracy. For example, the 
hole's mass and spin should be determined to < 0.1% ((J. It 
should even be possible to "map" the black hole's spacetime, 
testing whether it satisfies GR's stringent requirements 

Thanks to their extremal mass ratio, a formal prescription 
for modeling such systems now exists. At lowest order, the 
small body follows a geodesic orbit of the black hole [e.g., 
Ref. H, Eqs. (33.32a)-(33.32d)]. This must be corrected by 



the small body's interaction with its own spacetime distortion. 
The electromagnetic manifestation of this self force was given 
by Dirac 1 10], yielding the Abraham-Lorentz-Dirac equation 
of motion. Equations for a body's curved spacetime self in- 
teraction were worked out by Mino, Sasaki, and Tanaka full 
and by Quinn and Wald 1 12]; see Ref. [13] for an overview. 
Many [ 14] researchers are now working to develop practical 
schemes to compute the self force in black hole spacetimes, 
which is quite a challenge. 

Adiabatic radiation reaction (ARR). Astrophysical extreme- 
mass-ratio binaries allow a significant simplification for most 
of the inspiral, up to the last few orbits before the final plunge 
and merger. In this regime, the system evolves adiabatically: 
the radiation reaction time T ra d is much larger than the or- 
bital time T or b: T or b/T ra d ~ fJ>/M <C 1. The inspiral can 
be approximated as a flow through a sequence of geodesic or- 
bits. To compute the leading-order, adiabatic waveforms, it 
is necessary only to know the time-averaged rates of change 
of the three constants of motion: the energy E, axial angular 
momentum L z , and Carter constant Q fl flol . 

This approximation can be understood by expanding the 
small body's self acceleration a. Define e = fi/M. We can 
then write a = e [ag iss + ag ons + e (af ss + a.f ns ) +...]. 
Terms labeled "diss" describe dissipative aspects of the self 
force; they drive the inspiral. Those labeled "cons" are con- 
servative, representing non-dissipative components that con- 
tribute to the inspiraling body's inertia. Dissipative terms ac- 
cumulate secularly; conservative pieces do not. Their observ- 
able impact can be seen in the following expression obtained 
from a two-time expansion 1151 1 Kill for the azimuthal orbital 
phase (and correspondingly, the GW phase of each har- 
monic component of the waveform): 

$(t) = - [$ (i,ei)+E$i(i,rf) + O(£ 2 )] . (1) 

£ 

The leading-order, adiabatic waveforms contain only the term 
$o an d om it the subleading term $i which contributes a phase 
correction of order unity over the entire inspiral. The leading- 
order term is determined by aQ lss . Because it does not accu- 
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mulate secularly, ag ons does not contribute at leading order; it 
(along with a dlss ) contributes to the subleading term 

Estimates based on post-Newtonian theory suggest that the 
leading-order, adiabatic waveforms will be sufficiently accu- 
rate to detect waves in LISA data, and to provide an initial 
estimate of source parameters 1 17]. The phase mismatch be- 
tween $o (i) and the true signal is compensated for by adjust- 
ments in model parameters, which introduces systematic error 
in the inferred parameters. Eliminating this systematic error 
will require "measurement templates" that accurately model 
$i(i) — a far more difficult task. 

In this paper, we describe recent progress in computing adi- 
abatic waveforms, using (i) an approach called "Poor man's 
radiation reaction" (PMRR) for special classes of orbits, and 
(ii) augmenting PMRR using recent results of Mino 1 1 911 for 
fully generic orbits. 

Poor man 's radiation reaction. PMRR only requires knowl- 
edge of the energy and angular momentum the binary radi- 
ates; imposing global conservation, we evolve those quanti- 
ties, approximating inspiral as a flow through a sequence of 
orbits. PMRR cannot rigorously evolve all the constants E, 
L z and Q describing black hole orbits, except in special cases. 
It is simple to evolve E and L z — one computes the rate at 
which E and L z are carried to infinity 1 20] and absorbed by 
the hole's event horizon l2lll . and imposes global conserva- 
tion iH: E mb + E™ d = 0, L° rb + Lf d = 0. Carter's 
constant Q cannot be so evolved, since there is no notion of 
Q-flux carried by radiation and no associated global conser- 
vation law. Orbits with zero eccentricity or inclination are 
sufficiently constrained that Q is fixed by E and L z lE3ll . 

The general case admits no such constraints. Nevertheless 
we can produce inspirals accurate enough for data-analysis 
algorithm development using simple, crude approximations 
based on limiting cases. For example, forcing a certain incli- 
nation angle to be constant determines the inspiral using only 
E and L z fluxes [ 24]. Such inspirals are very unlikely to be 
accurate enough for detection templates. 

The central engine of PMRR is a complex function V>4 rep- 
resenting the small body's perturbation to the black hole's cur- 
vature. This function is the Weyl (vacuum) curvature tensor 
Cabcd projected onto a set of null vectors that are very conve- 
nient for describing outgoing radiation. Far away, ip4 is simply 
related to the two GW polarizations: 

1 d 2 

i/j 4 (r oo) = -— (h+ - ih x ) . (2) 

From this one finds the outgoing fluxes of E and L z 1 20]. This 
tp4 also completely describes the radiation's interaction with 
the black hole, and thus encodes the rate at which the large 
black hole absorbs ffandX, lEllEsll . 

Teukolsky derived a "master equation" for black hole per- 
turbations [26], and showed it separates by expanding in 



Fourier modes and spheroidal harmonics: 

= t 7^4 / ^Vli|„(r)S tau (%'(^ u( ) . 

(r — lacosvr J ^— ' 

v lm 

(3) 

(The parameter a = \S\/M is the black hole spin per unit 
mass.) It is not necessary to expand in modes: one can leave 
the Teukolsky equation as a PDE coupled in t, r, and 8 (the 
(j) dependence trivially decouples), and evolve initial data for 
?/>4- This works so well modeling source-free radiation lE7ll 
that it has been argued time domain methods may replace fre- 
quency decomposition for most applications 1 28 ] . For point 
particle sources, frequency domain methods are currently con- 
siderably more accurate [29]. 

The radial functions Rimui (r) are obtained by solving 

A 2 (A-X„J' - V lmu {r)R lmui = -T lmU) {r) , (4) 

where prime denotes d/dr, A = r 2 — 2Mr + a 2 , Vi muJ is apo- 
tential [e.g., [30], Eq. (4.3)], and % mU ) is a source constructed 
from the small body's stress energy tensor. We build a Green's 
function Gi mu (r, r') from solutions to the homogeneous wave 
equation (setting % moJ = 0), using certain analytic transfor- 
mations 13111 that allow high numerical accuracy. We then find 
Rlmu(r) by integrating Gi muJ (r, r') over the source, 

Rimu{r) = - f dr'G lmu (r,r')T lmu) (r') . (5) 

Our particular interest is in the limits r — > 00 and r — > the 
event horizon; from the function in these limits, we extract 
the E and L z fluxes mode-by-mode. Using many modes, we 
assemble the adiabatic rates of change (E) and (L z ). 

For bound orbits, the source is nicely described by a har- 
monic expansion. The continuous frequency ui goes over to a 
discrete set ui m k n = mSl^ + fcOa + nCl r , where ^l<j>,e,r de- 
scribe motion in 0, 6, and r Our modes become 4 
index objects, (Imkn). There is no mode-mode coupling, so 
this problem is ideally suited to parallel computation: differ- 
ent modes can be sent to different processors with little com- 
munication cost. Figure[2illustrates parallelization. The top 
panel shows the average CPU time per mode (Imkn) as a 
function of number of processors N . We find nearly linear 
scaling in 1/N, demonstrating that we have little paralleliza- 
tion overhead. Typical CPU time per mode is quite short, so 
that - 10 4 -10 5 modes can be computed in a reasonable time. 

The bottom panel shows the convergence of energy flux 
from a strong field orbit. We write the orbit's radial motion 
as r = pM/(l + ecosip); for this plot, we put eccentric- 
ity e = 0.5 and semi-latus rectum p = 4. The inclination 
l = 45°; the black hole's spin is a = 0.9 A/. (Precise defini- 
tions of p, e, and t can be found in l33ll ; their key property for 
this paper is that they are simply related to the constants E, 
L z , and Q.) We plot 

l K N 
El = ^ ^ ^ Elmkn (6) 
m=— i fc=— if n=—N 
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where Ei m k n is the modal energy flux. The cutoff values N 
and K (which formally are oc), are set large enough that ne- 
glected terms contribute a fractional amount less than 10~ 4 to 
the sum; details of this truncation will be presented elsewhere 
J34I1 . We generic ally find that E\ falls exponentially with I. 
The orbital parameter's influence on the rate of this falloff, and 
on other convergence criteria, will be presented in 1 34]. Most 
importantly, the convergence of these fluxes is fairly quick for 
e ^ 0.7, a very astrophy sic ally important range |5]. 
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FIG. 1: Top: Average CPU time per mode versus number of pro- 
cessors TV for two strong field orbits [circular (e = 0) and eccentric 
(e = 0.5)]. Our scaling shows T cpu oc 1/7V, showing how well 
PMRR parallelizes. The typical CPU time per mode is fairly short, 
though eccentric orbits are about 5-10 times slower than circular. 
Bottom: Flux dEi /dt (defined in the text) versus / for an orbit with 
p = 4, e = 0.5, and 1 = 45°; the black hole's spin a = 0.9M. We 
generically find that the flux sums converge exponentially with I. 

Mino 's adiabatic self force and Q 's evolution. Recent work by 
Mino 1 19] makes it possible to compute the time average (Q) 
of Q, allowing us to compute 3>rj(£) m Eq. Q for generic or- 
bits. The regularized self force f b i nvolves an integral over the 
orbiting body's past worldline II ill 1211 . Mino shows that this 
integral can be replaced in the adiabatic limit by a relatively 
simple expression involving the difference between "retarded" 
and "advanced" forces: f b = (/ t b ct - f b dv )/2. The "retarded" 
force depends on events on the orbiting body's past lightcone; 
the "advanced" force depends on the future lightcone. Their 
difference removes the divergent point-particle self interac- 
tion; the remaining force drives the inspiral. This result is 
strikingly similar to Dirac's result 1 10], and indeed reproduces 
the rule posited (without proof) by Gal'tsov 1 35]. 

The Carter constant is given by Q = Q a bP a p b , where Q a & is 
a Killing tensor and p a is the small body's 4-momentum. Tak- 
ing a time derivative yields Q = 2Q a bp a p b = 2Q a \ 3 p a j b . We 
express the self force in terms of the radiative Green's func- 
tion using Mino's result, and use the expansion of the radia- 
tive Green's function in terms of modes 1 35 ] . Time averaging 



yields an expression of the form jlrjl 

III />•// (r)]- (7) 

r—oo,H Imkn 

The quantities Rimkn(r) come from Ri mL0 (r) by writing 
(r)8(uj - oj mkn ). The first sum in 
Eq. means that the expression is evaluated near the horizon 
and near infinity. The quantities Qi m kn [t) are computed via 
an integral similar to that in Eq. 0, but with the source Ti muJ 
replaced by a new source built from the Killing tensor Q a \, and 
the orbiting body's 4-velocity, and evaluated at ui = u> m kn- Fi- 
nally W is the Wronskian which is independent of r near the 
horizon and near infinity as both Ri m kn ( r ) an d Qimkn{r) sat- 
isfy the homogeneous Teukolsky equation at those locations. 
Details of this calculation will be presented elsewhere 11611 . 
Using this result it will be as straightforward to compute (Q) 
as it is currently to compute (E) and (L z ). 
Applications and future directions. Building inspiral wave- 
forms requires us to compute backreaction effects upon a 
dense "grid" of orbits in parameter space. Each orbit is rep- 
resented by a coordinate (p, e, t); backreaction gives the tan- 
gent (p, e, i) to the inspiral trajectory. We flow along this tan- 
gent vector, building the trajectory \p(t), e(t), t(t)]. Choos- 
ing initial conditions, we likewise build the inspiral worldline 
z Ht) = [t,r(t),6(t),<f>(t)]. To date^only circular (e = 0) 
13611 and equatorial (l = 0,180°) |37] inspirals have been 
computed; the general case is under development 1 34], though 
we can present waveform "snapshots" (Fig.|2li. 



generic (eccentricity = 1/2, inclination = 45 degrees) 




200 



40(1 



600 



NIK) 



1000 



§■ 



circular (eccentricity = 0, inclination = 45 degrees) 



200 400 600 800 1000 



§ 



equatorial (eccentricity = 1/2, inclination = degrees) 






200 400 600 800 1000 

t/M 



FIG. 2: Gravitational waveform for several strong field orbits; all 
are for a hole with spin a = 0.9M, and have p = 4. Top panel: 
eccentricity e = 0.5, inclination t = 45°. Middle: e = 0, i = 
45°. Bottom: e = 0.5, t = 0°. The orbital parameters, particularly 
eccentricity, richly influence the wave's harmonic structure. 

The cases we have studied in detail demonstrate the rich 
physics encoded by these events. Three particularly interest- 
ing features are: 

• Tidal coupling to the event horizon: The orbiting body raises 
a tidal bulge on the black hole. This bulge interacts with the 
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orbit, transferring the hole's spin to the orbit, just as planetary 
tides can transfer angular momentum to a satellite. If the hole 
rotates rapidly and the orbit has shallow inclination, the tide 
can significantly prolong the inspiral |36]. 

• Strong field precessions: The time to oscillate through r 
does not equal the time to move through 2ir radians of <f>. 
The mismatch between these timescales gives perihelion pre- 
cession, a classic GR test. With black holes, this effect can 
amount to thousands of radians per orbit — the small body 
"whirls" many times near the hole before "zooming" out to 
large radius. This "zoom-whirl" behavior leaves a distinctive 
stamp on the waveform l37ll . seen in the short, high-frequency 
segments in Fig- EJ 

• Spin-orbit coupling: The black hole's spin makes the space- 
time geometry oblate, and drags spacetime into co-rotation 
with it. This splits the <f> and 8 frequencies, introducing fur- 
ther modulations to the waveform (middle panel of Fig.|2ji. 

All three of these features richly influence the phasing of a 
binary's GWs, and so are readily discernible in a phase coher- 
ent measurement. It is through determination of these features 
and their time evolutions that the binary's parameters and its 
strong field spacetime structure can be determined so well. 

Future work will develop ARR and waveforms with a focus 
upon templates for future GW searches. One goal is to use 
spectral methods for solving many of the formalism's equa- 
tions. Such methods are typically exponentially convergent in 
the number of basis functions. Since ARR requires many mul- 
tipoles, each must be as accurate as possible. We are very en- 
couraged by the work of Fujita and Tagoshi 1 38], who demon- 
strate that a particular set of basis functions allows modal so- 
lutions with essentially double precision accuracy. 
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